Bringing entanglement to the high temperature limit 
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We show the existence of an entangled nonequilibrium state at very high temperatures when two 
linearly coupled harmonic oscillators are parametrically driven and dissipate into two independent 
heat baths. This result has a twofold meaning: first, it fundamentally shifts the classical-quantum 
border to temperatures as high as our experimental ability allows us, and second, it can help 
increase by at least one order of magnitude the temperature at which current experimental setups 
are operated. 
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Introduction. — Since the establishment of quantum 
theory in last century there has been a long evolution 
on our concept of what is quantum and to what ex- 
tent it is required to explain observations in nature. At 
the very beginning the reduction postulate was proposed, 
clearly separating between quantum microscopic entities 
and classical macroscopic measuring apparatuses. Since 
then macroscopic quantum phenomena such as supercon- 
ductivity and coherent superposition in Bose-Einstein 
condensates [l[ , together with interference fringes of very 
massive molecules Q have been observed. Recently a 
proposal to create superpositions of dielectric bodies, 
such as viruses up to micron size, inside a high finesse 
optical cavity has been given Q. Hence the border be- 
tween the classical and quantum worlds seems to be more 
diffuse and intriguing than we could have conceived one 
century ago. 

Neither the usual transition criterion of h/S c h — > 
(with S c h the characteristic action of the system) is to 
be trusted, since this limit could be not completely con- 
tinuous and strong deviations have been reported in the 
semiclassical regime In the more realistic situation 
when the system interacts with the surrounding environ- 
ment, dissipation restricts purely quantum phenomena 
to within the very low temperatures limit [BJ, 

k B T/huj < 1 , (1) 

where ftw denotes the typical energy scale of the system 
and ksT the thermal energy. Above this limit, quantum 
correlations are inaccessible behind a 'mask' of thermal 
fluctuations. 

As a consequence, observing quantum phenomena im- 
plies the need for a very delicate pre-cooling process. But, 
is there any alternative to cooling for being quantum? 
In the present Letter, we defy the above classicality cri- 
terion and report the existence of a nonequilibrium en- 



tangled steady state for coupled harmonic oscillators at 
high temperatures, obtained through parametric driving. 
This result is quite fundamental, meaning that we might 
expect entanglement in hot highly nonequilibrium situa- 
tions, as pointed out @ for biological systems. Further, 
it could lighten the burden on quantum experiments re- 
quiring delicate pre-cooling setups. We note that though 
quantum coherence can play a role in biological processes 
at ambient temperature p}, demonstration of entangle- 
ment would be a much more extreme phenomenon. 

The model and its solution. — In order to render our 
arguments more quantitative, we study the entangle- 
ment between two interacting identical harmonic oscil- 
lators. Though an idealization, it encompasses a rea- 
sonable description of a wide variety of objects in nature 
such as nanomechanical oscillators [8[ , optical |9| and mi- 
crowave cavities [Toj . and movable mirrors [ll| to cite 
some, through which we expect to give a character of 
universality to the concepts that we expose here. The 
Hamiltonian of the system, "H s reads 

^ = E(?- + ^ mw2 ^) +c(<)0lQ2 ' (2) 

a = l \ 171 J 

with to the mass of the oscillator, uj the frequency and 
c(t) is the coupling coefficient. In what follows we as- 
sume: 

c(t) . 

= c + ci cos(w d t) , (3) 

TO 

that is, we consider a time dependent interaction, which 
plays a fundamental role in the creation and survival of 
entanglement. 

In any realistic scenario the system is not completely 
isolated from the outside. The most rigorous way to in- 
clude dissipation is by means of the system-bath model 
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[12j . We couple the oscillators to two independent baths 
(see figure QJ,), 

n = n s+ 2 f p ^ + ^l U, k - *%V ,(4) 

where the baths are modeled by an infinite collection 
of harmonic oscillators [HI, Q|. This independence ac- 
counts for a bath having a characteristic correlation 
length which exceeds the distance between oscillators. 
The opposite case corresponds to a model with common 
bath [1 5j j and leads to conservation of quantum entan- 
glement at higher than hoj/ks temperatures. This can 
be shown to be spurious since having different oscilla- 
tor frequencies or bath couplings (and no driving) leads 
again to the diagram in figure [Tb. Therefore we place 
ourselves in the most pessimistic situation for studying 
entanglement 

The evolution for the density matrix of the two oscil- 
lators, p s , can be cast as, 

P. (*>,*) = f d 4 X i J(X f ,t;X i ,0) Ps (X i ,0), (5) 

with X = {Q h+ ,Qi,-,Q 2 ,+ ,Q2,-} and J(X f ,t;X h 0) 
being the influence functional which is given in terms of 
a path integral expression after tracing out the environ- 
mental degrees of freedom [l7j]. Usually, the analytical 
evaluation of J(Xf, t; Xi, 0), even for time in dep endent 
systems, is only possible in very few cases [3 EH- Here, 
we have been able to derive an exact analytic expression 
for J(Xf, t; Xi, 0) in terms of the odd and even solutions 
of the Mathieu oscillator (See Appendices iBl and ICj) . The 
environmental influence enters via the spectral density 
I (to) = J2j Cj / (2nijOjj)5(ui — ojj). Here, we assume for 
simplicity Ohmic noise I(ui) = mjco. It produces white 
noise in the classical limit [18| ■ With our analytical result 
we can study the central system in any regime: low or 
high temperature, strong or weak damping, deeply quan- 
tum or semiclassical energy scales, etc., thus avoiding any 
spurious approximative corrections or limitations. As a 
result, any system that can be considered as two har- 
monic oscillators with linear coupling can be ascribed 
exactly to our description. 

Entanglement computation. — Linearity of the total 
Hamiltonian ensures that the state is always Gaussian, 
and thus its entanglement properties are fully character- 
ized by the covariance matrix crjj = + — 
(&)(£?') with £ = (Q\, Q2, -Pi) Pz)- An exact measure of 
entanglement is known for Gaussian states, the Logarith- 
mic Negativity En, as explained in Appendix |A"1 It is 
computed from the covariance matrix, which can be cal- 
culated from the propagator J(Xf, i; Xi, 0) (see appendix 
IE]) . In what follows we will exclusively use this measure. 

Entanglement in the time independent case. — In con- 
tact with an environment, each particle is asymptotically 
forced into a thermal state with a temperature equal to 
that of the bath it is connected to. This state is reached 



independently on the initial condition of the oscillator, 
which in the case of no driving [c\ = in ^} leads to 
the entanglement characteristics shown in figure ((Ud). 
That is, any state will, after thcrmalization, fall into 
either the blue (entangled) part or the white (separa- 
ble) part, depending only on the ratio co/mui 2 and the 
bath's temperature [16[ [30(. The entanglement region is 
restricted to the so called quantum limit fku < kT, as ex- 
pected from intuition, above such a temperature each os- 
cillator has an independent description because the quan- 
tum state is separable [3l| . 




FIG. 1: Generation of an entangled nonequilibrium state with 
dissipative environments, a, The system is formed by two lin- 
early coupled oscillators, initially thermalised due to each of 
them being dissipatively coupled to an environment at tem- 
perature T. Driving sinusoidally the coupling leads to produc- 
tion of entanglement even at very high temperatures, b, En- 
tanglement phase diagram for the case without driving. The 
state thermalizes to a state with no entanglement unless the 
temperature is below the quantum limit k^T < hco. c, Wigner 
phase-space representation of the normal modes. They are 
squeezed along orthogonal directions, so the oscillators are 
entangled (En — 0.33). The parameters are ksT/fioj = 10, 
7 = O.OIoj, ci = 0.5mo; 2 , while the snapshot has been taken 
at time tot = 6. 

Entanglement creation by driving — We sketch here a 
simple idea of how to produce an entangled nonequi- 
librium state at high temperatures. It may provide a 
huge leap in experimental requirements , while in addi- 
tion it definitely removes temperature from the list of 
possible criteria for classicality, the latter being an im- 
portant theoretical topic. The normal mode transfor- 
mation for the oscillator Hamiltonian ^ reads H = 
E Q =± P a/ 2m + <4<9<*/2 where Q± = (Q 1 ± Q 2 )/V2 
(P± = (Pi ± P 2 )/v / 2) and u? ± = oj ± c(t)/m. In the 
continuous variable setting, it is known that the max- 
imally entangled state -a kind of reference state, com- 
parison with which provides a quantification scheme for 
entanglement- is the Einstein, Podolsky, Rosen wavefunc- 
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tion [19J. It is just the infinite squeezing limit of the 
two-mode squeezed vacuum state, in which the indeter- 
minacies of Q+ and P_ are under the standard quantum 
limit set by Heisenberg's principle, while Q_ and J+ are 
above it (such that AQ±/AP± = exp (=F2r)/w 2 , with r 
the so-called squeezing parameter). The opposite situa- 
tion is also valid. Thus generation of entanglement can 
be provided by squeezing of the normal modes, which 
in turn can be generated through parametric driving of 
their frequencies ■ Both a time dependence in u or c 
will do, however the latter is better because it naturally 
provides the correct combination of squeezing between 
± modes. At the same time, the environment will try 
to destroy quantum coherence through equilibration to 
the thermal state. Thus we have two competing effects, 
whose balance will determine whether the steady state 
is entangled or not. In figure Q] we provide an example 
of normal mode squeezing in presence of the bath above 
the typical quantum limit ([T]) hsT/furf = 10 > 1. 

In figure [2] we summarize our results. Indeed, we find 
sets of parameters where entanglement is present at tem- 
peratures beyond the quantum limit, notice that in both 
figures kftT > fouj. Starting with a thermal state at the 
bath's temperature, the system reaches after a certain 
time a nonequilibrium steady state with nonzero entan- 
glement. We have chosen rather conservative couplings 
to the baths, as we will explain later, and still very high 
temperatures, k^T ^> fiu;, can be reached. 

It is a remarkable fact that while the system is forced 
into a highly nonequilibrium state, a steady state of en- 
tanglement is reached which is independent on the initial 
state of the system. To show this effect we plot in figure 
|3] (see inset) the time evolution of entanglement when the 
system starts with a two mode squeezed state and squeez- 
ing parameters r = 0,0.5, 1, and compare it to the case 
of an initial thermal state with the same temperature as 
the bath. 

New 'phase diagram' for entanglement — Parametric 
driving yields a new asymptotic behaviour which defines 
a new 'phase diagram', now dependent on four parame- 
ters: driving amplitude and frequency, temperature and 
the coupling to the bath. The driving frequency is overall 
chosen to be uod = 2 x 0.998w, and we also set Co = 0. 
While the optimal squeezing generation is obtained with 
a uid dependent on ui and c\ , the latter number seems to 
produce results nearly as good for different parameters, 
so it will be used unless otherwise stated. In figure [3] 
we see the points which delimit the border between pres- 
cncc(lcft)/absence(right) of entanglement, which is linear 
in temperature and driving amplitude and, as expected, 
the more isolated and driven the system is (low 7 and 
high ci), the higher the temperature can be reached. In 
addition to the exact result, we have plotted a simple 
estimation of the border which we explain next. 

We already mentioned that the entanglement produc- 
tion in this system can be viewed as a competition be- 
tween the squeezing due to the driving and mixing be- 
cause of the environment. The rate of squeezing can 




FIG. 2: (Color online) a, Time evolution of entanglement 
under parametric driving of the coupling at different en- 
vironmental temperatures ksT/hw = 250(red), 300(green), 
350(blue), with a damping of 7 = 0.005ojo, driving amplitude 
of ci = 0.5mwo and driving frequency u) c i = 2 x 0.998uo- A 
steady entangled state is reached in a reasonable time with 
a significative amount of entanglement, b, Now the temper- 
ature is kept fix , ksT/hjj = 5, with the same parameters, 
while the damping parameter is varied: 7 = 0.005ojo(red), 
0.01a;o (green), 0.02cj (blue). 



be obtained from the solutions to the nondissipative 
driven problem. They have the Mathieu form x(t) = 
exp (ifiMt)<p{t), where 4>{t) is a periodic function. If the 
Mathieu characteristic exponent finj is real, they are sta- 
ble, otherwise they are divergent which implies produc- 
tion of squeezing at a rate |Im(/ijy/)l (for every damped 
solution there is a divergent one) [21J. The rate of de- 
coherence can be estimated from the diffusion coefficient 
Z)[H (see Appendix iDj) . yielding jD ~ 7 ksT/fiiJ when- 
ever ksT > fvjj. Thus by comparison of both rates we 
obtain the new condition under which entanglement is 
present: 

ksT < |Im(/i M )| 

hw ~ 7 

which is seen to be a rather impressive match to the ex- 
act evolution. The condition above should be compared 
with the standard condition ([l} . In a nutshell the driving 
brings in a new quantum limit. 

Some examples — We give next some actual exam- 
ples of experiments which could profit from our strategy. 
However an additional comment is in order: the fact that 
squeezing grows approximately as |Im(/Ltjv/)|t also means 
that the energy and delocalization in space are increasing 
exponentially in time. Thus checking consistency with 
experimental size and energy considerations is a must. 
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FIG. 3: 'Phase diagram' of entanglement in the presence of 
parametric driving. We compare the condition ^ [lines] 
with the exact time evolution [dots] for different bath cou- 
plings 7 = 0.005a; (blue triangles), 0.001a; (green circles ) 
and 0.0005a; (red squares). Inset: time evolution for differ- 
ent initial conditions, namely a two mode squeezed vaccuum 
state (dotted curves) with squeezing parameter r — (red), 
0.5 (blue), 1 (green), as compared to that of an initial ther- 
mal state (black). They all converge after some tens of pe- 
riods. The parameters here are 7 = 0.001a;, ci = 0.2mo; 2 , 
u d = 2 x 0.9998oj and k B T/hu) = 10. 



Take for example two Calcium ions, each confined in 
its own planar Penning traps [22j ]. A trap can be fab- 
ricated by nanolithography with a size of d ~ 0.12/im. 
If a voltage of V = 10V is applied, the motional fre- 
quency is u> ~ 21 GHz and thus we can interpret fig- 
ure [T] as the temperature in Kelvin. A wire mediated 
capacitive coupling between traps allows to reduce the 
effective distance between ions and makes the coupling 
increase up to a reasonable level c(t) = cq = 0.047mwQ- 
If the frequencies are driven instead of the coupling (i.e. 
Lu(t) = luq + loi sinwdt), and assuming 7 = 0.0005a;, we 
still manage to get entanglement up to ~ 50K, while the 
derealization of the oscillators is yet below the trap size, 



ensuring no confinement leakage. To reach room temper- 
ature a very strong coupling would be required indeed, 
but our method allows the experimentalist to avoid build- 
ing a sub-4K (liquid Helium) setup. We believe this to 
be a huge experimental step. 

Another example is microwave superconducting cavi- 
ties [23| • The coupling between two cavities can be mod- 
ulated placing a superconducting qubit between them. 
The effective hamiltonian governing the dynamics is (|2|) . 
The typical frequencies in these resonators are in the Gi- 
gaHcrz regime, operating usually in the milikclvin range. 
The decoherence in these systems is 7 = 10 _4 o;, or even 
less. However the coupling is weak, around 10MHz. In 
this case, due to the weak coupling, the parametric driv- 
ing would enhance the amount entanglement that could 
be measured by nowadays technology [24 1. 

Current experiments with nanomechanical resonators 
have these typical parameters: u = 2ttv = 2ir x 15MHz, 
ci 



m = 10- 17 kg, ci - 10- 1 mw , and a quality factor Q 



20000, which yields a damping 7 = 5 x 10 



-5, 



An 



entangled state can be observed at IK. If the frequency 
can be increased a factor 10, it might reach the entangled 
regime in presence of liquid Helium. 

In addition it is notable that the strong coupling 
regime has been reached between a massive mechanical 
microresonator and light Q. Furthermore, a proposal for 
paramctrically driving the coupling between a nanome- 
chanical resonator and a superconducting electrical res- 
onator has been given in [2f| . Thus we might well foresee 
that these advances could be used to measure entangle- 
ment in yet unsuspected temperature regimes in the near 
future, while eliminating the need for complex and costly 
setups to cool objects to the quantum regime. 
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Appendix A: Entanglement quantification 



Entanglement can be easily quantified for a bipartite system of continuous variables in a Gaussian state. The 
logarithmic negativity (27j gives a characterization of the amount of entanglement which can be distilled into singlets. 
In the case of Gaussian continuous variable states, only the covariancc matrix is needed. The covariance matrix a is 
defined as 

=&&+6£i>/2 -<&>&> (Al) 
with £i = Q±, Q2, P\,Pi- The logarithmic negativity is defined as 

1 4 

E N = --^log2[Min(l,2|^|)] (A2) 

4=1 

where \ are the symplectic eigenvalues of the covariance matrix. They are simply the normal eigenvalues of the matrix 
— iEcr, with £ the symplectic matrix 



l a 
-1 2 



(A3) 



and I2 is the 2x2 identity matrix. 

Whenever the logarithmic negativity of the system is zero, we have a separable state p s = YliPiPi ® Pz\ an d eacn 
oscillator can be described independently. In continuous variable systems, the amount of entanglement is unbounded 
from above, having as a limiting case the maximally entangled EPR wavefunction with En — > 00. 

Appendix B: Decoupling the total system in normal modes 

The Hamiltonian of the total system reads 

H s = p- + \rm J j 2 q 2 + ^ + \mu 2 q 2 +c{t)q 1 q 2 , (Bl) 
2m 2 2m 2 

N 1 1 / \ 2 N 1 1 / / \ 2 

H IB = > -Z P k + o TO fe w fc x k o + / TT~rP k+ o m fe W k\ x k ; 77 , (B2) 

where c(t) = mco + mci cos(u;d£). Introducing the normal modes coordinates x+ and x^ defined by 

Qi = -j= (x+ + x-) , pi = (p+ +p~), (B3) 

92 = 71 ~ X ^ ' P2 = ^i=(P+-P-), (B4) 



H s reads 



H S = |+ + imf7 2 + (i)4 + ^ + imf7 2 _(i)x 2 _, 
2m 2 2m 2 



where = ± c(t)/m and ifjs 

W 1 . 1 1 c 2 

firs = JIo p k + o mkU3 k x l ~ -i=c k x k {x + + X-) + -(x + +x-) 2 

2m k 2 ^2 2V2m k ujj: 



N 



fc=l 

W 1 1 1 /2 



El /2 1 / /2 /2 1 / / / \ I c 'fe / \2 /r>r\ 

TTl-P k + 7; m ku k x k - —=c k x k {x+-x-)+ — 71 {x+ -X-) . (B6) 
2m'fc 2 y/2 2\J2m' k uj' k 



k=l 

These coordinates introduce a cross-term x+X- 1 which cancels out if 



m k w\ m' k oj' 2 k 



(B7) 
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This requirement does not means that the oscillators in the baths arc identic but their modes distributions. In the 
continuous limit, it implies that the spectral distributions characterizing the baths, Ji(w) and J2( w )j are the same. 
In this case, 



vJ 1 a 1 



-j^c k x k + -^=c fc a; fcj x+ - y-j=c k x k - —j=c k x kj x- I • (B8) 

This expression suggests the introduction of new set of coordinates qfc and Ofc defined by 

c\k = - 1 / - {ckx k + c'kx'k) , Qfc = - 1 7= (cfca; fc - c'fea;' fe ) , (B9) 
AfcV2 A fc v2 



which can be inverted 



£fc = ^^(A fe q fc +A fe Q fc ), a;',, = -=— {X k q k - A fc Q fc ) , (BIO) 
v2cfc v2c' fc 

-i- (A fcPfc + A k y k ) , p' k = -^—(X k p k -AM. (Bll) 



Pa 



After substituting in Hjb and choosing nikC 2 , = m'kc' k to eliminate a term proportional to p , we have 

A 2 m^ 2 
— q fc - A fc q fca;+ + 2^Vfc + 



tt _ / A fc h 2 i m kUk X l „2 x . _ , A fc m 2 | TO fc a; fc A fc 1 ^2 » „, \ / Rr v> 

^ - ^\2^ Pfc + ^^ qfe " AfcqfcX++ 2^ % + ^^ Q ^ A ^-/- (B12) 

To obtain a more standard version of the Hamiltonian, we could redefine m k — > A 2 /c 2 mfc and io\ = c 2 ,w 2 ,/\\ and 
impose A 2 = A 2 or just by choosing A 2 , = c\ = A 2 , so 

N f -. 2 i 2 

^ = E{^ + ^^ ± W+ + i^ + ! ^^ ±c ^-}- (B13) 

It means that we can conserve a small arbitrariness in the phase of the coupling by choosing different signs by Afc and 
Afc. However, for convenience we choose '+' for both. 
In summary, we have 

H = — + -mnW, + — + - m n 2 x 2 

2m 2 + + 2m 2 



N ( i 2 2 



fc=i 



or 



2-i ( 1 2 1 

* r i 

fc=l k K 



2 N * 9 

, LnSJ T 1 m2 , n»fcW fcrt a 



+ £ + 2 - - + E i — ^ + -r^a* - c *^- r • 

It is quite trivial, but we have derived an effective microscopic description of our initial assumption: normal modes 
coupled to identic but independent baths. It worths to be mentioned that not only the baths have the same modes, 

2 /2 

— fc-a- = — , k , but also the coupling between the system and the bath is the same, Afc = +Cfc = Afc. 

In order to complete our program an important point is left, if we want that the propagating function factorize, 
J\x-\-, x-,x' + , x'_] = J[x+, x', ] J[x—, x'_], obtaining that each normal mode evolves actually in an independent way, we 
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have to verify that the product by pairs of the equilibrium density matrix of the baths modes remains uncorrelatcd 
in the new coordinates. In the current case, the transformation of coordinates reads 



x k = -j={q k +Qk) , x' k = -j= (q fc - £j fc ) , 

Pk = ^= (Pfe + Vk) , P k = ^= (p fc - ¥k) , 
then, the product of the equilibrium density matrix of the fc-th mode of each bath reads 



m k uj k 



Z k \2TThsmh(uj k ni3) 



1 ( m'fcw'fe 



Z' k V27rfisinh(w / fc ft/3) 

m k uj k 
Z k \2Trhsmh(u) k h(3) 

m k u k 
Z k V27rftsinh(w fc n/3) 



cxp 



exp 



exp 



exp 



m k uj k 



2Trhsmh(uj k h(3) 



(( x i,k + x l fc) cosh ft/3) - 2x itk x i > tk ) 



((x'lf. + x'l k ) cosh(u}' k h/3) - 2x'i l kx' i ^k) 



2wHsmh(oj' kfi/3) 
m k uj k 



2nhsmh(io k hf3) 

m k u k 
2irhsmh{ui k h(3) 



((q?,fc + 4,k) cosh(w fe ^) - 2<\i, k <\v,k) 
{{Q 2 ik + Hi k ) cosh(oj k hp) - 20,j,Q,-j. i 



(B16) 
(B17) 



(B18) 



To obtain this desired result, we had to impose m k = m' k and u k = uj' k . So, it reduces our baths to be equal 
in detail, we mean, oscillator by oscillator. Only at this point we can affirm that the normal modes will evolve 
independently. This result for the bath modes can be interpret in geometrical terms as follows: the ispotcntial lines of 
two uncoupled identic harmonic arc defined by circumferences, so they are invariant under any rotation, which imply 
that the dynamical quantities obey exactly the same motion equations. It is important to mention that the normal 
modes are coupled to the bath in different coordinates than the real modes, however the introduction of the normal 
modes for the bath leaves the Jacobian of the transforation equals to 1, so after the trace the will generate completely 
equivalent results. 



Appendix C: Propagating function for the density matrix 

In normal modes, the evolution of the density matrix is governed by, 



p{ x +,f,y+j, x -,f,y-j,t) = J &z+,idy+,i J ^-,i^y-,iJ{x+j,y+j, x -j,y-j,t\x +!i ,y +!i ,x-^,y-^Q) 

x p( x +,i,y+d, x -,i,y-,i,t), (Ci) 

where J(x + j, y+j, X—j, y~j, t\x + j 1 y+.i, y-.i, 0) is the propagator of the reduced density matrix, 



J{ x +j,y+j^ x -j^y-j^\ x +,i^y+s^ x -,i^y-s^) 



Vx+ / Vy + / Vx^ / Vy_ 



exp ■ j j s ^ x+,x -^ ~ \F\ x +,y+-, x —,y-]i 



(C2) 



where S'[x + ,a;_] is the classical action and J r [x + ,y + ,X-,y-] the influence functional. T>x denotes an infinite product 
of measures in configuration space and implies a path integration over the paths x+(t), y+(t), X-(t) and y~{t) 
with endpoints x + (0) = x +ti , y(0) = y +ti , X-(0) = x~ :l , y_(0) = x + (t) = x + j, y{t) = y + j, X-{t) = X-j and 
?/_ (t) — y~j. However, at this point we have decoupled our system and we are describing it by two different harmonic 
oscillators coupled to identical but independent baths. So, 

p( x +J, y+J: X -J>y-J: *) = / dx+jdy+^dx-jdy- J+ (x+j , y+j ,t\x +ii , y +ti , 0) J_ (x_ ,/ , y_ j , t\x- ti ,y- ti , 0) 

x P(x+,i,y+,i,x- t i,y- t i,0), (C3) 



with 



J±( x ±,f,y±,f,t\ x ±,i,y±,i,0) = J T)x ± J T) y±cxp ^(S±[x±] - S±[y±})^ J 7 [x±,y±] 



(C4) 



For the case of a bath modeled by harmonic oscillators [13[, the general result for J 7 [x+,y+] was derived by Caldcira 
and Leggett [3] and it reads 



F[x+,y+} = exp 



h 2 



hi) / 

Jo 



(x+,i + y+,i) / ds^(s)[x+(s) -y+{s)} + 



Jo Jo 



dwy(s - u)[x+(u) + y+(u)][x+(s) - y+(s)\ 



x exp 



^JdsJ du[x+(u) - y+(u)]K(u - s)[x+(s) - y+(s)} 
similar expressions stands for the IF[x—,y-] mode, K(s) denotes the noise kernel 

K(s) = 



duo , / uih 

■ coth | — — — | cos(ujs)1(uj), 



2k B T 



(C5) 



(C6) 



wherein k B denotes the Boltzmann constant and T the temperature of the bath. The friction kernel j(s) in terms of 
the spectral density reads 



7(s) 



du I{ui) 



cos(ws), inOhmiccase j(s) = 2jS(s). 



(C7) 



An identical expression stands for T[x-,y-\. Since path integrals in J are quadratic, they can be done exactly to 
yield 



J 



1 



- exp \\{S+[xi\ S + [yf] + S_ [x c l] S_ [y% \ ^[xf, ^[xfi.yf 



(C8) 



being N± a normalization factor determined by the normalization of the propagator. To simplify further expressions, 
let's us to introduce the center of mass and difference variables, i.e., 



satisfying 



Q± = x±-y±, Q± = ^(x±-y±), 



Q±(s) - 79±(s) + ti±{s;tp)q±(s) = -2qf t ±~/6(t - s), 



Q ± (s) + jQ ± (s) + n 2 ± ( S ; <p)Q±(a) = -2Q h±1 S(s). 
It is important to mention that solution to these equations will be valid only for s > and it reads pll ] 

q±(a) = Vi ! ±(t,s;(f)q i> ±+v 2 ,±{t,s;cp)qf > ±, 



(C9) 

(CIO) 
(Cll) 

(C12) 



Q±(s) = Ui t ±(t,s;tp)Q it ±+U2,±(t,s;ip)Qf t ±. (C13) 
Since baths are defined by the same spectral density, then note that 7 is the same for ± cases. So we have that 

J( x +j>y+j, x -j,y-j> y+,i> %-,i,y-,i> °) = 



N(t) 



x exp 

x exp 

x exp 

x exp 



m{[b 3 ,+ (t; tp)q+,i - b 4i+ {t; tp)q +<f ]Q +tf + [b lt+ (t; <p)q+,i - b 2 ,+ (t; ip)q+j]} 
fn{[h-{t\<p)q-,i - b4,-(t;<p)q-j]Q+.f + [61,- (*; - b 2 -(t; <p)q- ,/]} 

{a n ,+(t\ y)q 2 +,i + [a\2,+ {t; (f) + o 2 i,+ (i; (f)}q + ,iq+j + a 22 . + (t; <p)q%j} 



-r{on,-(*;^)9-,i + [ai2,- (*;<£>) +a 2 i,-(t;(p)}q-,iq-,f + a 22 _{t]Lp)q 2 _ f } 



where N(t) = N+(t)N_(t), 



(Hi,± = \ [ ds i / ds 2 v it± (t 7 si;(p)vj ! ±(t,s 2 ;(p)K(s 1 - s 2 ), 
z Jo Jo 



(C14) 



(C15) 
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and 

&i,± (*;<£>) = «i,±(*)0;v) + 7i fe 2,± = wi !± (t,t;^), (C16) 
&3,±(*;¥>) = «2,±(*)0;<p), o 4 ,± = U2,±(t,i;y), (C17) 

Using last definitions we can express N± as N± = b ^(t) ■ Next step is the derivation of the master equation. We 
based our calculation on the paper of Zerbe and Hangii [2lj where the authors derived the exact quantum master 
equation for a single driven harmonic oscillator. 



Appendix D: Quantum Master Equation (QME) 



ift 



Quantum master equation for the normal modes of the initial system reads 



\-—\ 


' d 2 


2m \ 


K dx 2 + 


'--{ 


t Ql 


2m \ 


K dx 2 _ 


ih-y 

— (X+ 


-v+) 


ifvy 

— {X ~ 


-v-) 


-[D +<xp (t,0) 
m 



d 2 

d 2 ' 

dy 2 _ 



^\{t; V ){x\ 

Y^(*;v)0 



-vl) 



.-it) 



p(x+,y+,x-,y-) 

jxt or- j - - | ' 

( 3 3 \ 

V+) ( 77", q — J P(x+,y+,x-,y-) +iD +! p P (t,0)(x1 - y 2 _)p(x + ,y + ,X-,y-) 



- t^-^ P{x+,y+,x-,y- 



— [D-, xp (t,0) + D-, 
m 



y dx dy~ J Pv E +i y +i x -> y -> "l" M r > u A a 

f L» + . p J(a; + - y + ) + 7^- J p{x+,y+,x~, 



- LD_ 0) (x 2 . -y 2 _)p(x + ,y + ,x^,y^) 



where 



2/-) 
2/-), 



= At, 



. O 2 ,±04,± &4,± 
7 a 22,± - a22.± + 2t ^~ 

02, ± y 02, ±0? 

,*p(*,0) 



-T — ai2,± 

03, ± 



=,»(*, 0)= 2 ^6 4 ,± + ^ 

£>±jw(*,0) 

For small values of ft, D± tPX (t,Q) and D± 4>p (t, 0) can be written as [21 



'4,± 04, ± . 

7 °12,±, 

>3,± 03,± 

02, ± 

— : — fl i2,±- 

»2,±03,± 



J 2,± 
022,± + -, 



where A = ft 2 /3 2 /24m. 



D ± ,^,0) = ^ + ^(O? b (t)- 7 2 ) ) 
„ , , 2m7 2 A 



(Dl) 

(D2) 
(D3) 

(D4) 
(D5) 



Appendix E: Mean values and variances 



(f(x ± )) = J dQ ft± f(Q f>± )p(Q fi± ,q ft± = 0,t) (El) 



The first moments read in terms of the initial values (x±(to = 0)) = ((cc±,o)) and, 



(x ± (t)) = [/ 3 ,±(t) - |/i,±(t)]((x ±l0 )) + ^/i,±(*)<(p±,o)> (E2) 
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(p±(*)> = m-{x ± (t)) 



m[/ 2l ±(t) - |/i,±(t)]((x ±i0 )) + ^A,±(*)((p±,o)) 



(E3) 
(E4) 



The evolution of (p± (t)) is discontinuous at to = 0, i.e., lim t _j. + (p± (t)) = (p±.o) — Twy(x±fi) /2 es in general not equal 
to d>±,o}- This instantaneous jump of (p±{t)} can be removed with an environmental cutoff w c , or a non-factorizing 
initial state [l8|. The variances are obtained accordingly. They are given by 



<rx±x±(t) = [f2,± - o-ft.i) 



;/i,± (/a,± - |/i,±) 



±p± 



yi ±a ° p±p± +^fl ± a lh± , (E5) 



/2,±A,± — ^ (ji,±fo,± + A,±/2,± — ^/l,±A,± 



(ji,±h,± + A,±/2,± - 7A,±/i,±) o'SiPi + — 2 /i.±^p ±p± -+ 



(/l,±/ 



l,±Oll,± + /l,±Ctl2,± 



) , (E6) 



:(*) = ™ 2 (A.± - <x ± + 2mA,± (A,± - |A,±) < P± + AV: 





p±p± 



2hm 



(fi t ±au,± + 2A,±ai2,± + a22 : ± 



(E7) 



where we omitted the arguments of the functions ciij.± and fi,± for better lucidity. Here we note two missprints in 
[2l| . one is the presence of a global factor | in the last term of a x±p± and the other is in the last term of cr x±p± , in 

pl| it reads 2fi,m ^2/^ ^an,± + A.± a i2.± + ^22, ±) ■ Due to the discontinuity at t = 0, variances at i = + jump to 



^iii (^o+ ) 
(*0+) 



>x±p± 



Jp±p± 



x±x± ' 
-7CT + (7° 



(*o+) = l 2(J x ± x ± - 2 7 cr ° ±P± + CT ° ±P± . 



(E8) 
(E9) 
(E10) 



where i + means limt — > 0^ 



